Several coupled oscillators

Old stuff, to be rewritten. See the source in anotacioneslatex.tex.

They are important physical system what can be modeled, by approximation in their lagrangian to second order Taylor expansion, by

$$ M\ddot{X}=-K X $$

where $M$ and $K$ are symmetric matrices ($M$ comes from kinetic energy in the lagrangian and $K$ is the Hessian matrix of a potential function $V$ also in the exact lagrangian).

To solve it, we construct the associated first order system:

$$ Y=\left(\begin{array}{c} X \\ \dot{X} \end{array}\right) $$

and then

$$ \dot{Y}=\left(\begin{array}{cc} 0 & I \\ -M^{-1} K &0 \end{array}\right) \cdot Y $$ $$ \dot{Y}=A\cdot Y $$

To find eigenvalues and eigenvectors in order to use equation \ref{soluciondefinitiva} we observe that

\begin{itemize}

\item $\lambda$ is eigenvalue of $A$ iff $\lambda^2$ is eigenvalue of $-M^{-1}K$. This is so because the block matrix determinant formula

$$ det\left(\begin{array}{cc} H & B \\ C & D \end{array}\right)=det(H-BD^{-1}C)\cdot det(D)$$

applied to $\lambda I-A$.

\item If $$\xi=\left( \begin{array}{c} \psi \\ \eta \end{array}\right)$$

is an eigenvector of $A$ associated to $\lambda$ then $\psi$ is an eigenvector of $-M^{-1}K$ associated to $\lambda^2$. This is so because

$$ A\xi=\left(\begin{array}{cc} 0 & I \\ -M^{-1} K &0 \end{array}\right) \cdot \left( \begin{array}{c} \psi \\ \eta \end{array}\right)=\left( \begin{array}{c} \eta \\ -M^{-1}K\psi \end{array}\right)=\left( \begin{array}{c} \lambda \psi \\ \lambda \eta \end{array}\right) $$

So,

$$ A^2\xi=\left(\begin{array}{cc} -M^{-1} K & 0 \\ 0 &-M^{-1} K \end{array}\right) \cdot \left( \begin{array}{c} \psi \\ \eta \end{array}\right)=\left( \begin{array}{c} -M^{-1} K \psi \\ -M^{-1} K \eta \end{array}\right)= \left( \begin{array}{c} \lambda^2\psi \\ \lambda^2 \eta \end{array}\right) $$

And therefore

$$ -M^{-1} K \psi =\lambda^2\psi $$

Reciprocally, if $(\delta, \psi)$ is a pair such that $- M ^ { - 1 } K \psi = \delta \psi$ (of course $\delta$ is real and negative) then is easy to check that we have two eigenvalues with two eigenvectors for

$$ \left(\begin{array}{cc} 0 & I \\ -M^{-1} K &0 \end{array}\right) $: the pair $(\sqrt{\delta},\left( \begin{array}{c} \psi\\ \sqrt{\delta} \psi \end{array}\right))$ and the pair $(-\sqrt{\delta},\left( \begin{array}{c} \psi\\ -\sqrt{\delta} \psi \end{array}\right)$$

\end{itemize}

All the eigenvalues of $-M^{-1} K$ are reals, since is a symmetric matrix and they always can be diagonalized. Moreover, they are all negative numbers, because this matrix reflect a stable equilibrium point, i.e., from it potential energy increases in whatever direction we take. So we conclude:

\begin{itemize}

\item All the eigenvalues are pure imaginary ones. There is no $e^{\alpha t}$ terms.

\item For every eigenvalue of $A$ it can be chosen an eigenvector

$$ \psi=\left( \begin{array}{c} \xi \\ \eta \end{array}\right) $$ such that $\xi$ is real (because is an eigenvector of $-M^{-1} K$ and their eigenvalues are all real). So $$ Im(\psi)=\left( \begin{array}{c} 0 \\ Im(\eta) \end{array}\right). $$

\end{itemize}

So for every eigenvalue-eigenvector of $-M^{-1}K$, $(\lambda_n,\xi)$, we have two solutions:

$$ \Psi^{n}=e^{i\omega_n t} \xi^n $$

$$ \tilde{\Psi}^{n}=e^{-i\omega_n t} \xi^n $$

where $\omega_n=\frac{\sqrt{\lambda_n}}{i} \in \RR$ (remember $\lambda_n$ is negative). It may seem that this both solutions have the same initial value, which is a contradiction, but they have different initial velocities!

Moreover, observe that$

(\tilde{\Psi}^{n})^*=\Psi^n

$, because we can choose $\xi^n$ to be real.

So if we look for \textbf{real solutions} we have that

$$ \Psi(t)=\sum_n C_n \Psi^n+\tilde{C}_n\tilde{\Psi}^n $$

$$ (\Psi(t))^*=\sum_n C_n^* (\Psi^n)^*+(\tilde{C}_n)^*(\tilde{\Psi}^n)^* $$

must be equal, and since the solutions form a vector space we conclude that $C_n^*=\tilde{C}_n$

and so

$$ \Psi(t)=\sum_n C_n \Psi^n+C_n^*({\Psi}^n)^*=\sum_n 2Re[C_n\Psi^n]=\sum_n A_n cos(\omega_n t+ \phi_n) \xi^n $$

where $A_n$ and $\phi_n$, determined by initial conditions, are reals and have absorbed the constant to leave $\xi^n$ as a real vector.

\bigbreak

These fundamental solutions are called \textbf{normal modes}. The constants $\omega_k$ are called \textbf{fundamental frecuencies.}

A chain of a huge amount of coupled oscillators

\begin{center}

\includegraphics[width=15cm]{imagenes/noscillators.png}

\end{center}

We are going to study $N$ masses connected to a string with tension. We leave open ends, for the moment.

Let $\Psi_j (t)$ be the displacement from the equilibrium of every mass $j$. Let the mass be $m$ and the distance between the mass be $d$. Sometimes we will take $x=jd$ and write

$$ \Psi(x,t)=\Psi_{x/d} (t) $$

Analyzing every mass individually we obtain

$$ \ddot{\Psi}_1=-\frac{2T}{md}\Psi_1+\frac{T}{md}\Psi_2 $$

$$ \ddot{\Psi}_2=\frac{T}{md}\Psi_1-\frac{2T}{md}\Psi_2+\frac{T}{md}\Psi_3 $$

$$ \cdots $$

In general:

$$ \ddot{\Psi}=\frac{1}{md}\left( \begin{array}{cccc} {- T} & {T} & {0} & {0} \\ {T} & {-2 T} & {T} & {0} \\ {...} & {...} & {...} & ...\\ {0} & {0} & {T} & {-T} \end{array}\right) \cdot \Psi $$

or in short

$$ \ddot{\Psi}=-M^{-1}K\cdot \Psi $$

Observe that the matrix $-M^{-1}K$ is symmetric and all their eigenvalues are negatives, since it comes from a stable equilibrium point (see section \ref{coupledoscillators}). Moreover, from this section and the previous ones we know that a basis for solutions of this differential equation is

$$ A cos(\omega t+ \phi) \psi $$

where $-\omega^2$ is an eigenvalue of $-M^{-1}K$,

$\psi$ is an eigenvector for this eigenvalue, and $A$ and $\phi$ are real constants determined by the initial conditions: the $A$'s by the initial positions and the $\phi$'s by the initial velocities.

\bigbreak

\textbf{When $N$ is really big} is very difficult to find the eigenvalues by linear algebra, so we can proceed in a different way.

Also, the vector $\psi$ is a very large vector, since the number of masses $N$ will tend to infinity, and we want to study the shape'' of their components. In fact, when $N$ goes to infinity $\psi$ will be a function of position, giving the displacement of the masses in the initial time.

As usual, we will follow with the complex solution for $\Psi$ and take the real part in the final step. We think in solutions like

$$ \Psi_j(t)=e^{i\omega t} \xi(j) $$

where $\xi$ is a complex eigenvector of infinite length (keep an eye: in section \ref{coupledoscillators} we saw that with finite masses we can force $\xi$ to be real, but with infinite masses, a priori, we cannot).

To find infinite length eigenvectors we use the trick of \ref{symmetrytrick}. Since this infinite chain of oscillating masses has translation symmetry, we check that

$$ -M^{-1}K S=S (-M^{-1}K) $$

where

$$ -M^{-1}K= \left(\begin{array}{cccccc}{\ddots} & {\vdots} & {\vdots} & {\vdots} & {\vdots} & {\ddots} \\ {\cdots} & {-2 B} & {C} & {0} & {0} & {\cdots} \\ {\cdots} & {C} & {-2 B} & {C} & {0} & {\cdots} \\ {\cdots} & {0} & {C} & {-2 B} & {C} & {\cdots} \\ {\cdots} & {0} & {0} & {C} & {-2 B} & {\cdots} \\ {\ddots} & {\vdots} & {\vdots} & {\vdots} & {\vdots} & {\ddots}\end{array}\right) $$

and

$$ S=\left(\begin{array}{cccccc}{\ddots} & {\vdots} & {\vdots} & {\vdots} & {\vdots} & {\ddots} \\ {\cdots} & {0} & {1} & {0} & {0} & {\cdots} \\ {\cdots} & {0} &0& {1} & {0} & {\cdots} \\ {\cdots} & {0} & {0} & 0 &{1} & {\cdots} \\ {\cdots} & {0} & {0} & {0} & 0 & {\cdots} \\ {\ddots} & {\vdots} & {\vdots} & {\vdots} & {\vdots} & {\ddots}\end{array}\right) $$

What are the eigenvectors for $S$? We can take any $\beta \in \RR$, and then produce an eigenvector $\xi$. Observe that

$$ S\xi (j)=\xi(j+1)=\beta \xi(j) $$

so if we take $\xi(0)=1$ we conclude $\xi(j)=\beta^j$.

Moreover, since $j\in \ZZ$ we conclude that $|\beta|=1$ and therefore must be $\alpha \in [-\pi, \pi]$ such that $\beta=e^{i\alpha}$. In conclussion, for every $\alpha \in [-\pi, \pi]$ we have the eigenpair $(e^{i\alpha},\xi^{\alpha})$ where $\xi^{\alpha}(j)=e^{i\alpha j}$, respect to $S$.

But, what is the eigenvalue respect to $-M^{-1}K$?

$$ (-M^{-1}K \xi^{\alpha})(j) =-2 B e^{i\alpha j}+C e^{i\alpha (j-1)}+C e^{i\alpha (j+1)}= (-2B+Ce^{-i\alpha}+Ce^{i\alpha})e^{i\alpha j} $$

$$ =(-2B+C(e^{-i\alpha}+e^{i\alpha}))e^{i\alpha j}=(2C\cos(\alpha)-2B)e^{i\alpha j} $$

So $\xi^{\alpha}$ is an eigenvector of $-M^{-1}K$ associated to the eigenvalue $2C\cos(\alpha)-2B$, which must be negative (probably $C=B$, and we would get that).

Now, think that the pair $(2C\cos(\alpha)-2B,\xi^{\alpha})$ produces two solutions for the associated first order system (and then truncated''):

$$ \Psi^{\alpha}(t)=e^{i\omega_{\alpha} t} \xi^{\alpha} $$

$$ \tilde{\Psi}^{\alpha}(t)=e^{-i\omega_{\alpha} t} \xi^{\alpha} $$

where $\omega_{\alpha}=\frac{\sqrt{2C\cos(\alpha)-2B}}{i} \in \RR$. In the particular case where $B=C=\frac{T}{md}$,

$$ \omega_{\alpha}=2\sqrt{\frac{T}{md}}|\sin{\frac{\alpha}{2}}| $$

which is known as \textit{dispersion relation}

\bigbreak

A few remarks:

\begin{itemize}

\item The general solution is a linear combination of $\Psi^{\alpha}$'s and $\tilde{\Psi}^{\alpha}$'s. But since $\alpha$ has no restriction, beyond $\alpha \in [-\pi, \pi]$, because we don't have boundary conditions, $\alpha$ varies continuously. So, instead of a linear combination we get

$$ \Psi=\int_{-\pi}^{\pi} C(\alpha) \Psi^{\alpha} +\tilde{C}(\alpha) \tilde{\Psi}^{\alpha}d \alpha $$

Observe that this has been developed with a trick: we forgot the initial and final rows of the matrix $-M^{-1}K$ since we are dealing with a big $N$. This is the reason why we are finding more than $N$ eigenvalues (infinite, actually): since we are no taking into account the initial and final restriction, it is totally as if we had infinite masses.

\item Let's look for the general \textbf{real} solution. First, observe that:

$$ (\Psi^{\alpha})^*=e^{-i\omega_{\alpha} t} (\xi^{ \alpha})^*=e^{-i\omega_{\alpha} t} \xi^{ -\alpha}=\tilde{\Psi}^{-\alpha} $$

$$ (\tilde{\Psi}^{\alpha})^*=e^{i\omega_{\alpha} t} (\xi^{ \alpha})^*=e^{i\omega_{\alpha} t} \xi^{ -\alpha}={\Psi}^{-\alpha} $$

Then

$$ \Psi^*=\int_{-\pi}^{\pi} C^*(\alpha) (\Psi^{\alpha})^* +\tilde{C}^*(\alpha) (\tilde{\Psi}^{\alpha})^*d \alpha= \int_{-\pi}^{\pi} C^*(\alpha) \tilde{\Psi}^{-\alpha} +\tilde{C}^*(\alpha) {\Psi}^{-\alpha}d \alpha= $$

$$ =\int_{-\pi}^{\pi} C^*(-\alpha) \tilde{\Psi}^{\alpha} +\tilde{C}^*(-\alpha) {\Psi}^{\alpha}d \alpha $$

Since it must be $\Psi=\Psi^*$, and so :

$$ C^*(-\alpha)=\tilde{C}(\alpha) $$

(I don't know how to prove this!!, but the idea is that $\Psi^{\alpha}$ and $\tilde{\Psi}^{\alpha}$ constitute a basis for the vector space of solutions) and therefore:

$$ \Psi=\int_{-\pi}^{\pi} C(\alpha) \Psi^{\alpha} +\tilde{C}(\alpha) \tilde{\Psi}^{\alpha}d \alpha=\int_0^{\pi} 2Re[C(\alpha) \Psi^{\alpha} +\tilde{C}(\alpha) \tilde{\Psi}^{\alpha}]d\alpha $$

\item The parameter $\alpha$ is the seed for the wave number. Let $d$ be the distance between the masses, we can write $x=d\cdot j$ to be the position in the horizontal direction. If we try to write everything int terms of $x$ instead $j$ is useful to choose $\kappa=\frac{\alpha}{d}$ for the trial solution for the eigenvector. The normal modes would be:

$$ \Psi^{\kappa}(x, t)=a \cdot cos(\omega_{\kappa} t+\kappa x+\phi) $$

\item Dispersion relation is a name for the functional relation between wave number $\kappa$ and frecuency $\omega_{\kappa}$. In this case is

$$ \omega_{\kappa}=2\sqrt{\frac{T}{md}}sin(\frac{\kappa d}{2 }) $$

\item The idea for choosing $\xi(j)=A e^{i \alpha j}$ can be seen from other point of view.

The eigenvectors satisfy the relation

\begin{equation}\label{eigenvectorrelation}

\end{equation}

where we forget the first and the last relation because we assume $N$ very very large

One can observe that if $d$ is very small and taking into account the name change $x=jd$, equation \ref{eigenvectorrelation} can be interpreted (approximately) as:

$$ - m\omega^2 \xi(x)=T\left[ \frac{\xi(x+d)-\xi(x)}{d}-\frac{\xi(x)-\xi(x-d)}{d}\right] \approx T [\xi'(x+d/2)]-\xi'(x-d/2)] $$

and so

$$ - m\omega^2 \xi(x)\approx T d \xi''(x) $$

But since $m=d\mu$, where $\mu$ is the density and is supposed to be constant

$$ - \frac{\mu \omega^2}{T} \xi(x)\approx \xi''(x) $$

We have got an equation quite similar but in the variable $x$ or $j$ if you prefer, so is natural to choice $\xi(j)=A e^{i \alpha j}$.

\end{itemize}

\section{Boundary conditions}

We have studied an infinite system. If we want to come back to a finite system we impose boundary conditions, so we reduce the number of fundamental solutions that are allowed.

Particular case: $N=4$. We have:

$$ \ddot{\Psi}=\frac{1}{md}\left( \begin{array}{cccc}{- T} & {T} & {0} & {0} \\ {T} & {-2 T} & {T} & {0} \\ {0} & {T} & {-2T} & T\\ {0} & {0} & {T} & {-T}\end{array}\right) \cdot \Psi $$

We solve using Mathematica and find, evidently, four eigenvalues:

$$ -\omega_1^2=\frac{(-2-\sqrt{2}) T}{d m},-\omega_2^2=-\frac{2 T}{d m},-\omega_3^2= \frac{(-2+\sqrt{2}) T}{d m}, -\omega_4^2=0 $$

and four eigenvectors

$$ \xi_1=(-1,1+\sqrt{2},-1-\sqrt{2}, 1), \xi_2=(1,-1,-1,1), $$

$$ \xi_3=(-1,1-\sqrt{2},-1+\sqrt{2}, 1), \xi_4=(1,1,1,1) $$

With pictures ($T=2 ; m=1 ; d=0.5$):

\begin{tabular}{cc}

\includegraphics[width=60mm]{imagenes/sol1.png} &

\includegraphics[width=60mm]{imagenes/sol2.png} \\

\includegraphics[width=60mm]{imagenes/sol3.png} &

\includegraphics[width=60mm]{imagenes/sol4.png} \\

\end{tabular}

Could we obtain the same with the big-N-technique? Let's see. We start with infinite vibrating masses with displacement $\Psi_j(t)$, and so infinite fundamental frequencies

$$ \omega_{\alpha}=2 \sqrt{\frac{T}{m d}} \sin (\alpha / 2) $$

one for every $\alpha \in [0,\pi]$.

If we want to restrict to this particular case, we have to impose some conditions: \textbf{boundary conditions}. For example, if we force $\Psi_0=\Psi_1$ and $\Psi_4=\Psi_5$ we are removing the effect of masses 0 and 5 over the masses ranging from 1 to 4 (as desired). Let's study this two conditions: we are going to impose $\Psi_0=\Psi_1$ and $\Psi_4=\Psi_5$ and watch what $\alpha$'s survive.

$$ e^{i\omega_{\alpha}t} [ C(\alpha)+C(-\alpha)]=e^{i\omega_{\alpha}t} [ C(\alpha) e^{i\alpha}+C(-\alpha) e^{-i\alpha}] $$

$$ C(\alpha)+C(-\alpha)= C(\alpha) e^{i\alpha}+C(-\alpha) e^{-i\alpha} $$

and together with

$$ e^{i\omega_{\alpha}t} [ C(\alpha) e^{i\alpha 4}+C(-\alpha) e^{-i\alpha 4}]=e^{i\omega_{\alpha}t} [ C(\alpha) e^{i\alpha 5}+C(-\alpha) e^{-i\alpha 5}] $$

we arrive to

$$ \alpha=\frac{\pi n}{4} $$

Let's study this values one by one:

\begin{enumerate}

\item n=0. Then $\alpha=0$, $\omega_{\alpha}^2=0$ and $\xi=(1,1,1,1)$.

Correspond to the case $\omega_{4}$ following Mathematica.

\item $n=1$. Then $\alpha=\frac{\pi}{4}$, and $\omega_{\alpha}^2=\frac{T}{md}(2-\sqrt{2})$. We are in case 3 from Mathematica. The eigenvector obtained is $\xi=(e^{i\frac{\pi}{4}},e^{i\frac{\pi}{2}},e^{i\frac{3\pi}{4}},e^{i\pi})$, which is complex and do not coincide with the one obtained by Mathematica. But since $-\alpha$ produces the same eigenvalue, we can mix the previous eigenvector with $(e^{-i\frac{\pi}{4}},e^{-i\frac{\pi}{2}},e^{-i\frac{3\pi}{4}},e^{-i\pi})$ to make new eigenvectors for this eigenvalue. In fact is a (long) computation to check that the linear combination of this complex vectors that produces the same that the real above eigenvector is

$$ \left\{\left(x \rightarrow-\frac{1}{2}+i\left(-\frac{1}{2}+\frac{1}{\sqrt{2}}\right), y \rightarrow-\frac{1}{2}+i\left(\frac{1}{2}-\frac{1}{\sqrt{2}}\right)\right\}\right. $$

I found it by using Mathematica.

\item $n=2$. Then $\alpha=\frac{\pi}{2}$, and $\omega_{\alpha}^2=2\frac{T}{md}$. This correspond with case $\omega_2$ above. The eigenvector is not the same, but the phenomenom is the same that for $n=1$: we can recover by mixing both eigenvectors corresponding to $\alpha$ and $-\alpha$.

\item $n=3$. Then $\alpha=\frac{3\pi}{4}$, and $\omega_{\alpha}^2=\frac{T}{md}(2+\sqrt{2})$. Idem.

\item $n=4$. Then $\alpha=\pi$. Condition $\Psi_0^{\pi}=\Psi_1^{\pi}$ implies $C(\pi)+C(-\pi)=0$. Moreover, since $e^{i\pi j}=e^{-i\pi j}$, we deduce that $C(\pi)\Psi^{\pi}+C(-\pi)\Psi^{-\pi}=0$, so we can ignore this normal mode.

\item $n=5$. Then $\alpha=\pi+\frac{\pi}{4}$. As we reasoned above, this would be the same case as $\pi+\frac{\pi}{4}-2\pi=-\frac{3 \pi}{4}$ which, in fact, is paired with $\alpha=\frac{3\pi}{4}$ (case $n=3$)

\item All the remaining cases are repeated

\end{enumerate}

________________________________________

________________________________________

________________________________________

Author of the notes: Antonio J. Pan-Collantes

antonio.pan@uca.es


INDEX: